/*
** ShpTETR4.c  -  this file contains methods for the 4-noded linear
**                tetrahedron solid finite element shape.
**
**    The adopted numbering scheeme is shown below:
**
**      The nodes are numbered according to the following rule:
**      1.  first, counterclockwise, the nodes at the face with 
**          natural coordinate L4 = 0,
**      2.  and finally the node at the vertice with L4 = 1.
**
**                            + 4
**                           /|\
**                          / | \
**                         /  |  \
**                        /   |   \
**                       /    |    \
**                      /     |     \
**                     /      |      \
**                    /       |       \
**                   /        |        \
**                  /         |         \
**                 /          |          \
**               1+ . . . . . | . . . . . +2
**                 \_         |         _/
**                   \_       |       _/
**                     \_     |     _/
**                       \_   |   _/
**                         \_ | _/
**                           \|/
**                            + 3
**
**
**           Internal numbering:
**
**
**                            3  (L4=1)
**                           /|\
**                          / | \
**                         /  |  \
**                        /   |   \
**                       /    |    \
**                      /     |     \
**                     /      |      \
**                    /       |       \
**                   /        |        \
**                  /         |         \
**                 /          |          \
**        (L1=1)  0 . . . . . | . . . . . 1  (L2=1)
**                 \_         |         _/
**                   \_       |       _/
**                     \_     |     _/
**                       \_   |   _/
**                         \_ | _/
**                           \|/
**                            2  (L3=1)
**
** --------------------------------------------------------------------------
**
** Entry Points:
**
** --------------------------------------------------------------------------
**
** void InitTETR4( void )
**
**      Initializes sub-class TETR4.
**
** --------------------------------------------------------------------------
**
** static Telement_3d *NewTETR4( void )
**
**      shp           -  pointer to element shape              ( out )
**
**      This method allocates memory for a given TETR4 element shape
**      and returns a pointer to the new object.
**
** --------------------------------------------------------------------------
**
** static void ReadReseqTETR4 ( int *reseq )
**
**      reseq    -  node numbering resequecing vector          (  in )
**
**      Reads node numbering resequecing vector TETR4 finite
**      element shape class.
**
** --------------------------------------------------------------------------
**
** static void ReseqTETR4 ( int *reseq )
**
**      reseq    -  node numbering resequecing vector          (  in )
**
**      Specifies node numbering resequecing for TETR4 finite
**      element shape class.
**
** --------------------------------------------------------------------------
**
** static void GetReseqTETR4 ( int *reseq )
**
**      reseq    -  node numbering resequecing vector          ( out )
**
**      Gets resequence vector for TETR4 finite element shape node
**      numbering.
**
** --------------------------------------------------------------------------
**
** static void ConnectTETR4 ( Telement_3d *shp,
**                            int    *connect )
**
**      shp      -  given element shape                        (  in )
**      connect  -  connectivity of the given element shape    ( out )
**
**      Returns connectivity of a given TETR4 finite element shape.
**
** --------------------------------------------------------------------------
**
** static void DerivMapTETR4rst ( Telement_3d    *shp,
**                                NatCoord  *gauss_rst,
**                                DerivNat  deriv_rst[4] )
**
**      shp           -   given element shape  (not used)      (  in )
**      gauss_rst     -   natural coords.                      (  in )
**      deriv_rst     -   partial derivative wrt r,s,t         ( out )
**
**      This routine  returns the  "map"  functions  derivatives  with
**      respect to  the natural coordinates for a 4-noded  tetrahedron
**      solid finite element shape. The "map" function derivatives are
**      evaluated at the given natural coordinates. Numbering scheme is 
**      same as that given above.
**
** --------------------------------------------------------------------------
**
** static int  DimTETR4 ( Telement_3d *shp )
**
**      shp          -  given element shape  (not used)        (  in )
**
**      Returns dimension number of TETR4 shape, ie. 3.
**
** --------------------------------------------------------------------------
**
** static void GaussTypeTETR4 ( Telement_3d    *shp,
**                              GaussType *gauss_type )
**
**      shp        -  given element shape (not used)             (  in )
**      GaussType  -  returned Gauss integration type            ( out )
**
**      Returns the Gauss integration type which is consistent with
**      this element shape, i.e., TETR_QUADRATURE.
**
** --------------------------------------------------------------------------
**
** static int  GetFirstEdgeTETR4 ( Telement_3d    *shp,
**                                 int       *num_verts,
**                                 int       connect[] )
**
**      shp        -  given element shape                        (  in )
**      num_verts  -  number edge verts                          ( out )
**      connect    -  edge connectivity                          ( out )
**
**      Gets no. of verts and connectivity of first edge of given shape.
**      The verts are returned in sequential order along the edge.
**
** --------------------------------------------------------------------------
**
** static int  GetNextEdgeTETR4 ( Telement_3d    *shp,
**                                int       *num_verts,
**                                int       connect[] )
**
**      shp        -  given element shape                        (  in )
**      num_verts  -  number edge verts                          ( out )
**      connect    -  edge connectivity                          ( out )
**
**      Gets no. of verts and connectivity of next edge of given shape.
**      The verts are returned in sequential order along the edge.
**      If there is no other edge returns a false (0) status.
**      Otherwise, it returns a true (1) status.
**
** --------------------------------------------------------------------------
**
** static int  GetFirstFaceTETR4 ( Telement_3d    *shp,
**                                 int       *num_verts,
**                                 int       connect[] )
**
**      shp        -  given element shape                        (  in )
**      num_verts  -  number face verts                          ( out )
**      connect    -  face connectivity                          ( out )
**
**      Gets no. of verts and connectivity of first face of given shape.
**      The verts are returned in sequential order along the face.
**
** --------------------------------------------------------------------------
**
** static int  GetNextFaceTETR4 ( Telement_3d    *shp,
**                                int       *num_verts,
**                                int       connect[] )
**
**      shp        -  given element shape                        (  in )
**      num_verts  -  number face verts                          ( out )
**      connect    -  face connectivity                          ( out )
**
**      Gets no. of verts and connectivity of next face of given shape.
**      The verts are returned in sequential order along the face.
**      If there is no other face returns a false (0) status.
**      Otherwise, it returns a true (1) status.
**
** --------------------------------------------------------------------------
**
** static int  GetShpEdgeTETR4 ( Telement_3d    *shp,
**                               int       corners[],
**                               ShapeType *type,
**                               int       *num_nodes,
**                               int       connect[] )
**
**      shp        -  given element shape                        (  in )
**      corners    -  corners of edge                            (  in )
**      type       -  returned edge shape type                   ( out )
**      num_nodes  -  number edge nodes                          ( out )
**      connect    -  edge connectivity                          ( out )
**
**      Gets shape type, no. of nodes, and connectivity of an edge of
**      given shape for a given pair of corner nodes.
**      The returned edge shape type is LINE2.
**
** --------------------------------------------------------------------------
**
** static int  GetShpFaceTETR4 ( Telement_3d    *shp,
**                               int       corners[],
**                               ShapeType *type,
**                               int       *num_nodes,
**                               int       connect[] )
**
**      shp        -  given element shape                        (  in )
**      corners    -  corners of face                            (  in )
**      type       -  returned face shape type                   ( out )
**      num_nodes  -  number face nodes                          ( out )
**      connect    -  face connectivity                          ( out )
**
**      Gets shape type, no. of nodes, and connectivity of a face of
**      given shape for a given triple of corner nodes.
**      The returned face shape type is TRIA3.
**
** --------------------------------------------------------------------------
**
** static void GetShpDomainTETR4  ( Telement_3d    *shp,
**                                  ShapeType *type,
**                                  int       *num_nodes,
**                                  int       connect[] )
**
**      shp        -  given element shape                        (  in )
**      type       -  returned domain shape type                 ( out )
**      num_nodes  -  number domain nodes                        ( out )
**      connect    -  domain connectivity                        ( out )
**
**      Gets shape type, no. of nodes, and connectivity of given shape.
**      The returned shape type is TETR4  and the number of domain nodes
**      is 4.
**
** --------------------------------------------------------------------------
**
** static void MakeTETR4 ( Telement_3d *shp,
**                         int    *connect )
**
**      shp      -  given element shape                        (  in )
**      connect  -  connectivity of the given element          (  in )
**
**      Makes a TETR4 element shape of given connectivity.
**
** --------------------------------------------------------------------------
**
** static void MakeDefaultTETR4 ( Telement_3d *shp,
**                                int    *connect )
**
**      shp      -  given element shape                        (  in )
**      connect  -  connectivity of the given element          (  in )
**
**      Makes a TETR4 element shape of given connectivity.
**      The connectivity is assumed in the default class sequence.
**
** --------------------------------------------------------------------------
**
** static void MapFuncTETR4 ( Telement_3d    *shp,
**                            NatCoord  *nat_coord, 
**                            ShapeVal  *shape )
**
**      shp          -  given element shape (not used)         (  in )
**      nat_coord   -   natural coords.                        (  in )
**      shape       -   shape funcs. evaluated at given pt.    ( out )
**
**      This routine returns the "map" functions  for a  4-noded solid
**      tetrahedron  element shape. The "map" functions are  evaluated
**      at the given natural coordinates. 
**
** --------------------------------------------------------------------------
**
** static void NodalNatCoordTETR4 ( Telement_3d    *shp,
**                                  NatCoord  *coord )
**
**      shp         -   given element shape                    (  in )
**      coords      -   natural coordinates of its nodes       ( out )
**
**      Returns shape natural nodal coordinates.
**
** --------------------------------------------------------------------------
**
** static int  NumNodesTETR4 ( Telement_3d *shp )
**
**      shp          -  given element shape  (not used)        (  in )
**
**      Returns number of nodes in a 4-noded tetrahedron shape, ie. 4.
**
** --------------------------------------------------------------------------
**
** static int  ReadTETR4 ( Telement_3d *shp, int *conn )
**
**      shp           -  Given element shape                   (  in )
**      conn          -  Connectivity of element               (  in )
**      
**      This method reads nodal connectivity  data of a given "TETR4" 
**      shape element.
**
** --------------------------------------------------------------------------
**
** version: 1-000
**
** history: Created:   30-Jul-96   Luiz Fernando Martha
**              Stolen from Femoop and modified.
**
**          Modified:  19-Oct-98   Joaquim Bento Cavalcante Neto
**              Modified some methods to agree with other shape modules.
**
** Modified:  Novem/2000 - Antonio Miranda (amiranda)
**      Changed the input of read method (FILE -> int *conn)
**
** --------------------------------------------------------------------------
*/

/* Global variables and symbols           */

#define _SHPTETR4_C

#include <stdio.h>
#include <stdlib.h>

#include "shp.h"


/*  Local variables and symbols  */

/*
** --------------------------------------------------------------------------
** Local data structures:
**
*/

typedef struct _datatetr4 {
   int   node[4];    /* element node indices */
} DataTETR4;

/* Resequence node numbering vector:
 */
static int  rs[4] = {  0,  1,  2,  3  };

/* Use quadrature 4 gauss point */
static int num_gauss_point = 4;
static NatCoord gauss_rst[4] = {
     { 0.5854101966249685, 0.1381966011250150, 0.1381966011250150 },
     { 0.1381966011250150, 0.5854101966249685, 0.1381966011250150 },
     { 0.1381966011250150, 0.1381966011250150, 0.5854101966249685 },
     { 0.1381966011250150, 0.1381966011250150, 0.1381966011250150 }
};
static GaussWgt gauss_wgt[4] = { 1.0/24.0, 1.0/24.0, 1.0/24.0, 1.0/24.0 };

static ShapeVal shape_rst[4][4];
static DerivNat deriv_rst[4][4];


/* Edge and Face lists:
 */
#define NEDGES  6
static int curr_edge = 0;;
static int el[NEDGES][2] =  { {   0,  1  },
                              {   1,  2  },
                              {   2,  3  },
                              {   0,  2  },
                              {   0,  3  },
            {   1,  3  } };

#define NFACES  4
static int curr_face = 0;;
static int fl[NFACES][3] =  { {   0,  1,  2  },
                              {   0,  2,  3  },
            {   0,  3,  1  },
            {   2,  1,  3  } };

/* Local entry points: */

static Telement_3d *NewTETR4        ( void );
static void ReadReseqTETR4     ( int * );
static void ReseqTETR4         ( int * );
static void GetReseqTETR4      ( int * );
static void ConnectTETR4       ( Telement_3d *, int * );
static void DerivMapTETR4rst   ( Telement_3d *, NatCoord *, DerivNat[]  );
static int  DimTETR4           ( Telement_3d * );
static void MapFuncTETR4_i     ( Telement_3d *, int , ShapeVal ** );

static int    NumGaussPoint      (void);
static void   DerivMapTETR4_i   ( Telement_3d *, int, DerivNat **);
static double GaussWgtTETR4_i   ( Telement_3d *shp, int i);

static int  GetFirstEdgeTETR4  ( Telement_3d *, int *, int [] );
static int  GetNextEdgeTETR4   ( Telement_3d *, int *, int [] );
static int  GetFirstFaceTETR4  ( Telement_3d *, int *, int [] );
static int  GetNextFaceTETR4   ( Telement_3d *, int *, int [] );
static int  GetShpEdgeTETR4    ( Telement_3d *, int [], AuxShapeType *, int *, int [] );
static int  GetShpFaceTETR4    ( Telement_3d *, int [], AuxShapeType *, int *, int [] );
static void GetShpDomainTETR4  ( Telement_3d *, ShapeType *, int *, int[] );
static void MakeTETR4          ( Telement_3d *, int * );
static void MakeDefaultTETR4   ( Telement_3d *, int * );
static void MapFuncTETR4       ( Telement_3d *, NatCoord *, ShapeVal * );
static void NodalNatCoordTETR4 ( Telement_3d *, NatCoord * );
static int  NumNodesTETR4      ( Telement_3d * );
static int  ReadTETR4          ( Telement_3d *, int *conn );

/* Public Entry Points: */


/* ============================= InitTETR4 ================================== */

void InitTETR4( void )
{

   shape_class[NF_TETR4].New                   =    NewTETR4;
   shape_class[NF_TETR4].readresequence        =    ReadReseqTETR4;
   shape_class[NF_TETR4].resequence            =    ReseqTETR4;
   shape_class[NF_TETR4].getresequence         =    GetReseqTETR4;
   shape_class[NF_TETR4].connectivity          =    ConnectTETR4;
   shape_class[NF_TETR4].shapefunc             =    MapFuncTETR4;
   shape_class[NF_TETR4].shapefunc_i           =    MapFuncTETR4_i;

   shape_class[NF_TETR4].derivmaprst           =    DerivMapTETR4rst;
   shape_class[NF_TETR4].derivshprst           =    DerivMapTETR4rst;
   shape_class[NF_TETR4].derivmap_i            =    DerivMapTETR4_i;
   shape_class[NF_TETR4].derivshp_i            =    DerivMapTETR4_i;

   shape_class[NF_TETR4].derivxyz              =    shp_derivXYZ;
   shape_class[NF_TETR4].numgausspoint         =    NumGaussPoint;
   shape_class[NF_TETR4].bmatrix               =    shp_BMatrix;
   shape_class[NF_TETR4].gaussw_i              =    GaussWgtTETR4_i;

   shape_class[NF_TETR4].dim                   =    DimTETR4;
   shape_class[NF_TETR4].free                  =    shp_freeshape;
   shape_class[NF_TETR4].gausstonode           =    shp_gausstonode;
   shape_class[NF_TETR4].jacobian              =    shp_jacobian;

   shape_class[NF_TETR4].getfirstedge          =    GetFirstEdgeTETR4;
   shape_class[NF_TETR4].getnextedge           =    GetNextEdgeTETR4;
   shape_class[NF_TETR4].getfirstface          =    GetFirstFaceTETR4;
   shape_class[NF_TETR4].getnextface           =    GetNextFaceTETR4;
   shape_class[NF_TETR4].getshpedge            =    GetShpEdgeTETR4;
   shape_class[NF_TETR4].getshpface            =    GetShpFaceTETR4;
   shape_class[NF_TETR4].getshpdomain          =    GetShpDomainTETR4;
   shape_class[NF_TETR4].getverts              =    ConnectTETR4;
   shape_class[NF_TETR4].isinterface           =    shp_isinterface;
   shape_class[NF_TETR4].localsys              =    shp_localsyssolid;
   shape_class[NF_TETR4].make                  =    MakeTETR4;
   shape_class[NF_TETR4].makedefault           =    MakeDefaultTETR4;
   shape_class[NF_TETR4].mapfunc               =    MapFuncTETR4;
   shape_class[NF_TETR4].nodalnatcoord         =    NodalNatCoordTETR4;
   shape_class[NF_TETR4].nodetovertex          =    shp_nodetovertex;
   shape_class[NF_TETR4].nummapnodes           =    NumNodesTETR4;
   shape_class[NF_TETR4].numshpnodes           =    NumNodesTETR4;
   shape_class[NF_TETR4].numverts              =    NumNodesTETR4;
   shape_class[NF_TETR4].print                 =    shp_printshape;
   shape_class[NF_TETR4].read                  =    ReadTETR4;

  
   int i;
   for (i = 0; i < num_gauss_point; i++)
   {
     MapFuncTETR4 (NULL, &gauss_rst[i], shape_rst[i]);
     DerivMapTETR4rst ( NULL, &gauss_rst[i], deriv_rst[i]);
   }   
}


/* ============================== NewTETR4 ================================== */

static Telement_3d *NewTETR4 ( void )
{
   Telement_3d      *shp;

   shp = (Telement_3d *)calloc( 1, sizeof(Telement_3d) );

   shp->type = NF_TETR4;

   return( shp );
}


/* ============================ ReadReseqTETR4 ============================== */

static void ReadReseqTETR4 ( int *reseq )
{
   int i;

   for (i=0;i<4;i++)
    rs[i] = reseq[i] - 1;
}


/* ============================== ReseqTETR4 ================================ */

static void ReseqTETR4 ( int *reseq )
{
   int i;

   for( i=0; i<4; i++ )
    rs[i] = reseq[i];
}


/* ============================= GetReseqTETR4 ============================== */

static void GetReseqTETR4 ( int *reseq )
{
   int i;

   for( i=0; i<4; i++ )
     reseq[i] = rs[i];
}


/* =============================== ConnectTETR4 ============================= */

static void ConnectTETR4 ( Telement_3d *shp,
                           int    *connect )
{
   connect[0] = shp->inc[0];
   connect[1] = shp->inc[1];
   connect[2] = shp->inc[2];
   connect[3] = shp->inc[3];
}


/* ========================== DerivMapTETR4rst ============================= */

static void DerivMapTETR4rst ( Telement_3d   *shp,
                               NatCoord *gauss_rst,
                               DerivNat deriv_rst[4] )
{
   double  L1, L2, L3, L4;

   ignore (shp);

   L1 = gauss_rst->r;
   L2 = gauss_rst->s;
   L3 = gauss_rst->t;
   L4 = 1.0 - L1 - L2 - L3;

   deriv_rst[0].r =  1.0;
   deriv_rst[1].r =  0.0;
   deriv_rst[2].r =  0.0;
   deriv_rst[3].r = -1.0;

   deriv_rst[0].s =  0.0;
   deriv_rst[1].s =  1.0;
   deriv_rst[2].s =  0.0;
   deriv_rst[3].s = -1.0;

   deriv_rst[0].t =  0.0;
   deriv_rst[1].t =  0.0;
   deriv_rst[2].t =  1.0;
   deriv_rst[3].t = -1.0;
}


/* =============================== DimTETR4 ================================= */

static int DimTETR4( Telement_3d *shp )
{ 
   ignore (shp);

   return ( 3 );
}


/* ============================== NumGaussPoint =========================== */
static int NumGaussPoint (void)
{
  return (num_gauss_point);
}

/* ============================== DerivMapTETRA4_i =========================== */
static void DerivMapTETR4_i  ( Telement_3d *shp, int i, DerivNat **ptr)
{
  if (i < 0 || i > 3)
  {
    ptr = NULL;
    return;
  }
  *ptr = deriv_rst[i];
}


/* ============================== GaussWgtTETRA4_i =========================== */
static double GaussWgtTETR4_i  ( Telement_3d *shp, int i)
{
  if (i < 0 || i > 3)
    return 0.0;

  return (gauss_wgt[i]);
}



#if 0
/* ============================== GaussTypeTETR4 ============================ */

static void GaussTypeTETR4                            GaussType *gauss_type )
{
   ignore (shp);

   *gauss_type = TETR_QUADRATURE;
}

#endif


/* ============================= GetFirstEdgeTETR4 ========================= */

static int  GetFirstEdgeTETR4  ( Telement_3d    *shp,
         int       *num_verts,
                int       connect[] )
{
   *num_verts = 2;
   curr_edge = 0;

   connect[0] = shp->inc[el[curr_edge][0]];
   connect[1] = shp->inc[el[curr_edge][1]];
   return( 1 );
}


/* ============================= GetNextEdgeTETR4 =========================== */

static int  GetNextEdgeTETR4 ( Telement_3d    *shp,
             int       *num_verts,
             int       connect[] )
{

   *num_verts = 2;
   curr_edge++;

   if( curr_edge >= NEDGES )
    return( 0 );

   connect[0] = shp->inc[el[curr_edge][0]];
   connect[1] = shp->inc[el[curr_edge][1]];
   return( 1 );
}


/* ============================= GetFirstFaceTETR4 ======================== */

static int  GetFirstFaceTETR4 ( Telement_3d    *shp,
        int       *num_verts,
        int       connect[] )
{

   *num_verts = 3;
   curr_face = 0;

   connect[0] = shp->inc[fl[curr_face][0]];
   connect[1] = shp->inc[fl[curr_face][1]];
   connect[2] = shp->inc[fl[curr_face][2]];
   return( 1 );
}


/* ============================= GetNextFaceTETR4 =========================== */

static int  GetNextFaceTETR4  ( Telement_3d    *shp,
        int       *num_verts,
        int       connect[] )
{
   *num_verts = 3;
   curr_face++;

   if( curr_face >= NFACES )
    return( 0 );

   connect[0] = shp->inc[fl[curr_face][0]];
   connect[1] = shp->inc[fl[curr_face][1]];
   connect[2] = shp->inc[fl[curr_face][2]];
   return( 1 );
}


/* =============================== GetShpEdgeTETR4 ========================== */

static int  GetShpEdgeTETR4 ( Telement_3d    *shp,
                              int       corners[],
                              AuxShapeType *type,
                              int       *num_nodes,
                              int       connect[] )
{
   int          status = 0;

   *type = LINE2;
   *num_nodes = 2;

   if( corners[0] == shp->inc[0] )
   {
    connect[0] = shp->inc[0];
    if( corners[1] == shp->inc[1] )
    {
     connect[1] = shp->inc[1];
     status = 1;
    }
    else if( corners[1] == shp->inc[2] )
    {
     connect[1] = shp->inc[2];
     status = 1;
    }
    else if( corners[1] == shp->inc[3] )
    {
     connect[1] = shp->inc[3];
     status = 1;
    }
   }
   else if( corners[0] == shp->inc[1] )
   {
    connect[0] = shp->inc[1];
    if( corners[1] == shp->inc[0] )
    {
     connect[1] = shp->inc[0];
     status = 1;
    }
    else if( corners[1] == shp->inc[2] )
    {
     connect[1] = shp->inc[2];
     status = 1;
    }
    else if( corners[1] == shp->inc[3] )
    {
     connect[1] = shp->inc[3];
     status = 1;
    }
   }
   else if( corners[0] == shp->inc[2] )
   {
    connect[0] = shp->inc[2];
    if( corners[1] == shp->inc[0] )
    {
     connect[1] = shp->inc[0];
     status = 1;
    }
    else if( corners[1] == shp->inc[1] )
    {
     connect[1] = shp->inc[1];
     status = 1;
    }
    else if( corners[1] == shp->inc[3] )
    {
     connect[1] = shp->inc[3];
     status = 1;
    }
   }
   else if( corners[0] == shp->inc[3] )
   {
    connect[0] = shp->inc[3];
    if( corners[1] == shp->inc[0] )
    {
     connect[1] = shp->inc[0];
     status = 1;
    }
    else if( corners[1] == shp->inc[1] )
    {
     connect[1] = shp->inc[1];
     status = 1;
    }
    else if( corners[1] == shp->inc[2] )
    {
     connect[1] = shp->inc[2];
     status = 1;
    }
   }

   return( status );
}


/* =============================== GetShpFaceTETR4 ========================== */

static int  GetShpFaceTETR4 ( Telement_3d    *shp,
                              int       corners[],
                              AuxShapeType *type,
                              int       *num_nodes,
                              int       connect[] )
{
   int          status = 0;

   *type = TRIA3;
   *num_nodes = 3;

   if( corners[0] == shp->inc[0] )
   {
    connect[0] = shp->inc[0];
    if( (corners[1] == shp->inc[1]) && (corners[2] == shp->inc[2]) )
    {
     connect[1] = shp->inc[1];
     connect[2] = shp->inc[2];
     status = 1;
    }
    else if( (corners[1] == shp->inc[3]) && (corners[2] == shp->inc[1]) )
    {
     connect[1] = shp->inc[3];
     connect[2] = shp->inc[1];
     status = 1;
    }
    else if( (corners[1] == shp->inc[2]) && (corners[2] == shp->inc[3]) )
    {
     connect[1] = shp->inc[2];
     connect[2] = shp->inc[3];
     status = 1;
    }
   }
   else if( corners[0] == shp->inc[1] )
   {
    connect[0] = shp->inc[1];
    if( (corners[1] == shp->inc[2]) && (corners[2] == shp->inc[0]) )
    {
     connect[1] = shp->inc[2];
     connect[2] = shp->inc[0];
     status = 1;
    }
    else if( (corners[1] == shp->inc[0]) && (corners[2] == shp->inc[3]) )
    {
     connect[1] = shp->inc[0];
     connect[2] = shp->inc[3];
     status = 1;
    }
    else if( (corners[1] == shp->inc[3]) && (corners[2] == shp->inc[2]) )
    {
     connect[1] = shp->inc[3];
     connect[2] = shp->inc[2];
     status = 1;
    }
   }
   else if( corners[0] == shp->inc[2] )
   {
    connect[0] = shp->inc[2];
    if( (corners[1] == shp->inc[0]) && (corners[2] == shp->inc[1]) )
    {
     connect[1] = shp->inc[0];
     connect[2] = shp->inc[1];
     status = 1;
    }
    else if( (corners[1] == shp->inc[1]) && (corners[2] == shp->inc[3]) )
    {
     connect[1] = shp->inc[1];
     connect[2] = shp->inc[3];
     status = 1;
    }
    else if( (corners[1] == shp->inc[3]) && (corners[2] == shp->inc[0]) )
    {
     connect[1] = shp->inc[3];
     connect[2] = shp->inc[0];
     status = 1;
    }
   }
   else if( corners[0] == shp->inc[3] )
   {
    connect[0] = shp->inc[3];
    if( (corners[1] == shp->inc[1]) && (corners[2] == shp->inc[0]) )
    {
     connect[1] = shp->inc[1];
     connect[2] = shp->inc[0];
     status = 1;
    }
    else if( (corners[1] == shp->inc[2]) && (corners[2] == shp->inc[1]) )
    {
     connect[1] = shp->inc[2];
     connect[2] = shp->inc[1];
     status = 1;
    }
    else if( (corners[1] == shp->inc[0]) && (corners[2] == shp->inc[2]) )
    {
     connect[1] = shp->inc[0];
     connect[2] = shp->inc[2];
     status = 1;
    }
   }

   return( status );
}


/* ============================== GetShpDomainTETR4 ========================= */

static void GetShpDomainTETR4 ( Telement_3d    *shp,
                                ShapeType *type,
                                int       *num_nodes,
                                int       connect[] )
{
   *type = NF_TETR4;
   *num_nodes = 4;

   ConnectTETR4( shp, connect );
}


/* ================================ MakeTETR4 ============================== */

static void MakeTETR4 ( Telement_3d *shp,
                        int    *connect )
{
   shp->inc[rs[0]] = connect[0];
   shp->inc[rs[1]] = connect[1];
   shp->inc[rs[2]] = connect[2];
   shp->inc[rs[3]] = connect[3];
}


/* ================================ MakeDefaultTETR4 ======================= */

static void MakeDefaultTETR4 ( Telement_3d *shp,
                               int    *connect )
{

   shp->inc[0] = connect[0];
   shp->inc[1] = connect[1];
   shp->inc[2] = connect[2];
   shp->inc[3] = connect[3];
}


/* ============================= MapFuncTETR4 =============================== */

static void MapFuncTETR4 ( Telement_3d   *shp,
                           NatCoord *nat_coord, 
                           ShapeVal *shape )
{
   double  L1, L2, L3, L4;

   ignore (shp);

   L1 = nat_coord->r;
   L2 = nat_coord->s;
   L3 = nat_coord->t;
   L4 = 1.0 - L1 - L2 - L3;

   shape[0] = L1;
   shape[1] = L2;
   shape[2] = L3;
   shape[3] = L4;
}

/* ============================= MapFuncTETR4_i =============================== */

static void MapFuncTETR4_i (Telement_3d *shp, int i, ShapeVal **shape )
{
  if (i < 0 || i > num_gauss_point)
  {
    *shape = NULL;
    return;
  }
  *shape = shape_rst[i];
}

/* =========================== NodalNatCoordTETR4 =========================== */

static void NodalNatCoordTETR4 ( Telement_3d    *shp,
                                 NatCoord  *coord )
{
   ignore (shp);

   coord[0].r =  1.0; coord[0].s =  0.0; coord[0].t =  0.0;
   coord[1].r =  0.0; coord[1].s =  1.0; coord[1].t =  0.0;
   coord[2].r =  0.0; coord[2].s =  0.0; coord[2].t =  1.0;
   coord[3].r =  0.0; coord[3].s =  0.0; coord[3].t =  0.0;
}

/* =========================== NodalCoordTETR4 ============================== */

static void NodalCoordTETR4 ( Telement_3d    *shp,
                              NodeCoord *coord )
{
#if 0
   int      i;
   DataTETR4   *dat;

   dat = (DataTETR4 *) shp->dat;

   for (i=0;i<4;i++) {
      coord[i].x = node_coords[shp->inc[i]-1].x;
      coord[i].y = node_coords[shp->inc[i]-1].y;
      coord[i].z = node_coords[shp->inc[i]-1].z;
   }
#endif
}


/* ============================ NumNodesTETR4 =============================== */

static int NumNodesTETR4( Telement_3d *shp )
{
   ignore (shp);

   return ( 4 );
}


/* ============================= ReadTETR4 ================================== */

static int ReadTETR4( Telement_3d *shp, int *conn )
{
   int i;

   for (i=0;i<4;i++)
     shp->inc[rs[i]] = conn[i];

   return( TRUE );
}


/* =========================================================== End of File */
